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Abstract 

Liapunov  functions  are  constructed  for  nonlinear  systems  of  ordinary  differential  equa¬ 
tions  whose  linearized  system  at  an  equilibrium  point  possesses  cither  a  simple  zero  eigen¬ 
value  or  a  complex  conjugate  pair  of  simple,  pure  imaginary  eigenvalues.  The  construction 
is  explicit,  and  yields  parametrized  families  of  Liapunov  functions  for  such  systems.  In  the 
case  of  a  zero  eigenvalue,  the  Liapunov  functions  contain  quadratic  and  cubic  terms  in  the 
state.  Quartic  terms  appear  as  well  for  the  case  of  a  pair  of  pure  imaginary  eigenvalues. 
Predictions  of  local  asymptotic  stability  using  these  Liapunov  functions  are  shown  to  co¬ 
incide  with  those  of  pertinent  bifurcation- theoretic  calculations.  The  development  of  the 
paper  is  carried  out  using  elementary  properties  of  multilinear  functions.  The  Liapunov 
function  families  thus  obtained  are  amenable  to  symbolic  computer  coding. 
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1  INTRODUCTION 


In  this  paper,  we  construct  families  of  Liapunov  functions  useful  in  assessing  the  asymptotic  stability 
of  critical  equilibrium  points  of  a  class  of  systems 

i  =  /(*)  (1) 

where  x  £  IRC  and  /  is  at  least  four  times  continuously  differentiable.  Throughout  the  presentation, 
we  let  the  origin  x  —  0  be  the  equilibrium  point  of  interest  of  this  system  of  ordinary  differential 
equations.  Critical  cases  in  the  study  of  stability  of  the  origin  of  Eq.  (1)  are  those  in  which  the 
Jacobian  matrix  Df[ 0)  possesses  at  least  one  eigenvalue  with  zero  real  part,  and  no  eigenvalues 
with  positive  real  part.  In  these  situations,  it  is  not  possible  to  ascertain  whether  or  not  the  origin 
of  (1)  is  locally  asymptotically  stable  solely  from  the  linearized  system  x  =  (D/( 0))x. 

We  focus  on  two  specific  critical  cases  in  the  stability  analysis  of  Eq.  (1).  These  correspond 
to  the  Jacobian  matrix  of  (1)  at  the  origin  possessing  either  a  simple  zero  eigenvalue  or  a  complex 
conjugate  pair  of  simple,  pure  imaginary  eigenvalues.  The  associated  hypotheses  ((S)  and  (II) 
below,  respectively)  are  introduced  next.  The  techniques  of  this  paper  can  also  be  applied  to  other 
critical  cases,  such  as  those  involving  multiple  critical  eigenvalues  [7]. 

In  essence,  the  paper  relies  on  two  basic  tools  in  achieving  its  goals.  First,  the  notation  of 
multilinear  functions  is  adhered  to  throughout  the  paper,  in  denoting  terms  in  the  Taylor  series 
expansions  of  the  nonlinear  system  of  interest  as  well  as  Liapunov  function  candidates  and  their 
derivatives.  Second,  a  new  result  on  local  definiteness  of  a  class  of  scalar  bivariate  functions 
is  introduced.  This  appears  as  Lemma  1  in  Section  4.  Lemma  1  is  indeed  key  to  the  subsequent 
constructions  of  Liapunov  functions  in  Cases  (S)  and  (II).  These  constructions  are  rather  mechanical 
given  Lemma  1,  albeit  somewhat  tedious. 

The  first  critical  case  of  interest  in  this  work  is  characterized  by  the  occurrence  of  a  simple 
zero  eigenvalue  of  Df( 0),  with  all  remaining  eigenvalues  having  strictly  negative  real  parts.  As 
discussed  in  [2],  stability  of  the  origin  in  this  situation  is  closely  related  to  the  stability  of  bifurcated 
equilibrium  points  in  smooth  parametrizations  of  Eq.  (1).  Because  of  this  connection  to  stationary 
(or  static)  bifurcation,  this  critical  case  will  be  referred  to  here  as  “Case  (S).” 

( S )  The  Jacobian  Df( 0)  possesses  a  simple  zero  eigenvalue,  with  all  other  eigenvalues  in  the  open 
left  half  of  the  complex  plane. 
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In  the  second  critical  case  of  interest  here,  D/( 0)  is  assumed  to  possess  a  single  complex  con¬ 
jugate  pair  of  simple,  pure  imaginary  eigenvalues,  with  the  remaining  eigenvalues  lying  in  the  open 
left  half  complex  plane.  Reference  [1]  discusses  the  relationship  of  stability  of  the  origin  in  this  case 
with  stability  of  bifurcated  periodic  solutions  of  smooth  parametrizations  of  (1).  Hopf  bifurcation 
to  periodic  solutions  occurs  for  generic  such  parametrizations  of  Eq.  (1)  under  these  circumstances. 
Because  of  this  connection  to  Hopf  bifurcation,  this  critical  case  will  be  referred  to  as  “Case  (H).” 

(H)  The  Jacobian  0/(0)  possesses  a  complex  conjugate  pair  of  simple,  pure  imaginary  eigenvalues, 
with  all  other  eigenvalues  in  the  open  left  half  of  the  complex  plane. 

Recently,  there  has  been  significant  interest  in  feedback  stabilization  of  nonlinear  systems  in 
critical  cases  (see,  e.g.,  the  review  paper  [18]  and  references  therein).  This  has  yielded  various 
existence  and  synthesis  results  on  stabilizability  by  either  smooth  or  continuous  feedback.  The  main 
contribution  of  the  present  paper  to  this  body  of  work  is  the  construction  of  new  Liapunov  functions 
for  such  systems,  in  the  two  critical  cases  (S)  and  (H).  Liapunov  functions  facilitate  estimation  of 
the  domain  of  attraction  of  a  stable  equilibrium  point,  and  as  such  can  serve  to  quantify  the 
efficacy  of  a  given  control  design.  Liapunov  functions  can  be  used  to  define  performance  indices 
in  optimization-based  feedback  control  design  of  nonlinear  systems.  Such  performance  indices  can 
involve  estimates  of  the  achieved  domain  of  attraction  and  measures  of  the  adequacy  of  the  transient 
response.  The  Liapunov  functions  derived  here  for  Cases  (S)  and  (H)  are  given  explicitly  in  terms 
of  the  system’s  dynamics,  and  are  amenable  to  symbolic  computer  coding.  As  a  by-product  of  our 
results,  known  formulae  for  testing  stability  in  the  critical  cases  (S)  and  (II)  (so-called  bifurcation 
formulae')  are  found  to  follow  easily  from  the  Liapunov  functions  we  obtain.  Note,  however,  that 
these  formulae  alone  do  not  yield  analytical  performance  indices  of  the  type  just  alluded  to. 

The  asymptotic  stability  of  nonlinear  systems  in  critical  cases  has  received  significant  attention 
in  the  literature  (e.g.,  [3],  [11],  [13]-[15],  [17],  [19]).  Although  some  of  these  works  have  employed 
Liapunov  stability  analysis,  the  Liapunov  functions  used  have  generally  been  defined  only  implicitly. 
In  some  cases,  this  is  linked  to  the  use  of  implicitly  defined  nonlinear  coordinate  transformations 
to  lower  dimensional  problems.  Implicitly  defined  Liapunov  functions  suffice  when  the  goal  of  the 
analysis  is  limited  to  deriving  sufficient  conditions  for  local  asymptotic  stability.  For  instance, 
one  result  of  Mees  and  Chua  [17]  gives  a  Liapunov  function  for  planar  systems  (1)  satisfying 
(H).  Implications  of  this  result  for  higher  dimensional  models  (1)  follow  from  the  Center  Manifold 
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Theorem  (cf.  [17], [5]).  In  the  present  paper,  we  give  an  explicit  construction  of  families  of  Liapunov 
functions  for  critical  nonlinear  systems  satisfying  either  hypothesis  (S)  or  (H)  which  apply  directly 
to  the  given  n-dimensional  system  description  (1). 

This  paper  is  organized  as  follows.  In  Section  2,  pertinent  results  on  multilinear  functions 
are  given.  The  set-up  for  construction  of  Liapunov  functions  for  systems  (1)  is  formulated  using 
multilinear  function  notation  in  Section  3.  A  lemma  giving  sufficient  conditions  for  local  definiteness 
of  a  class  of  scalar  bivariate  functions  is  presented  in  Section  4.  Section  5  contains  a  result  on 
solutions  of  Liapunov  matrix  equations  for  a  coefficient  matrix  with  a  zero  eigenvalue  or  a  pair  of 
pure  imaginary  eigenvalues.  The  main  results  of  the  paper  appear  in  Sections  6  and  7.  Section 
6  contains  an  explicit  construction  of  a  family  of  Liapunov  functions  for  Case  S,  and  Section  7 
contains  an  analogous  construction  for  Case  (H).  Conclusions  are  collected  in  Section  8. 

Notation.  In  what  follows,  IRn  denotes  the  space  of  n-dimensional  column  vectors  having  real 
entries,  while  (Pn  denotes  the  space  of  n-dimensional  column  vectors  with  complex  entries.  The 
complex  conjugate  of  a  quantity  (scalar,  vector,  or  matrix)  a  is  denoted  by  a.  The  transpose  of  a 
vector  or  matrix  a  is  denoted  aT .  For  a  vector  space  V,  denote  by  ( V)k  the  vector  space  obtained  as 
the  fc-tuple  product  V  x  ••  •  x  V.  The  Jacobian  derivative  of  a  function  4>  is  denoted  D<j>.  The  norm 
of  a  vector  x  £  IRn  will  be  denoted  |aj,  and  the  same  notation  will  apply  to  any  compatible  matrix 
norm.  Denote  by  r  (resp.  /)  the  right  column  (resp.  left  rowr)  eigenvector  of  Df( 0)  corresponding 
to  the  critical  eigenvalue  0  (Case  (S))  or  iuic  (Case  (H)).  For  consistency  with  previous  literature  [1], 
[2],  [6],  [9],  the  first  component  of  r  is  set  to  unity,  and  /  is  then  chosen  subject  to  the  normalization 
Ir  =  1.  (Ensuring  that  the  first  component  of  r  is  nonzero  in  some  cases  requires  a  reordering  of 
the  elements  of  x.) 
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2  RESULTS  ON  MULTILINEAR  FUNCTIONS 


Multivariable  Taylor  series  can  be  conveniently  represented  in  terms  of  multilinear  functions.  We 
shall  employ  multilinear  functions  in  representing  Taylor  series  expansions  both  for  the  vector  field 
/(x)  of  Eq.  (1),  and  for  the  Liapunov  functions  whose  construction  is  the  main  purpose  of  this 
work.  In  this  section,  we  present  several  useful  facts  pertaining  to  multilinear  functions. 

2.1  Multilinear  Functions 

Multilinear  functions  may  be  defined  as  follows. 

Definition  1.  Let  Vi,  V2,. . .  ,14  and  W  be  vector  spaces  over  the  same  field.  A  map  ip  :  V\  x  I2  x 
.  . .  x  14  t— >  W  is  said  to  be  multilinear  (or  k-hnear )  if  it  is  linear  in  each  of  its  variables.  That  is 
[4,  p.  76],  for  arbitrary  v'  ,vl  £  Vt,  i  —  1, . . . ,  fc,  and  for  arbitrary  scalars  a,  a,  we  have 

^(i>V  •  • ,  avl  +  av\  . . .  ,vk)  =  aij>(vl . . ,  v* , . . . ,  vk)  -f  aip(vl , . . . ,  v1 , . . . ,  vk).  (2) 

□ 

We  refer  to  k  as  the  degree  of  the  multilinear  function  ip.  In  particular,  multilinear  functions  of 
degree  two,  three  and  four  are  referred  to  as  bilinear,  trihnear  and  teirahnear  functions,  respectively. 

We  shall  in  the  sequel  deal  exclusively  with  multilinear  functions  ip  whose  domain  is  the  product 
space  of  k  identical  vector  spaces  V\  —  14  =  . . .  =  14  =  V.  For  such  multilinear  functions,  we  have 
the  following  notion  of  symmetry. 

Definition  2.  A  fc-linear  function  ip  :  V  x  V  x  . . .  x  V  m  IE  is  symmetric  if,  for  any  vl  £  V , 
i  =  l, ... ,k,  the  vector 

ip{vl ,  v2,...,vk)  (3) 

is  invariant  under  arbitrary  permutations  of  the  argument  vectors  v1.  □ 

With  an  arbitrary  multilinear  function  ip,  we  associate  a  symmetric  multilinear  function  i/’s 
resulting  from  the  following  simple  device,  known  as  the  symmetrization  operation  [4,  pp.  88-89]. 
Given  a  multilinear  function  ip(x*  ,x2,..  .  ,xk ),  define  a  new  (symmetric)  multilinear  function  i/'s  as 
follows: 

ips(x:  ,x2, . . .  ,xk)  :=  —  ip{x11  ,x'2  , . . .  ,xlk),  (4) 

where  the  sum  is  taken  over  the  id  permutations  of  the  integers  1,2,  ...,k. 
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2.3  Coordinate  Representation  of  Scalar  Multilinear  Functions 

In  this  subsection,  we  state  a  useful  representation  result  for  scalar  multilinear  functions  x p  : 
(lRn)^  t— >  lRm .  The  representation  rests  upon  a  choice  of  basis  (“coordinates”)  for  IRn .  Thus, 
let  {?-1 ,  r2, . . . ,  r"}  be  a  basis  for  IRn .  By  a  standard  result  [4,  Proposition  3.6.1],  to  this  basis  there 
corresponds  a  unique  dual  basis  which  we  may  view  as  consisting  of  row  vectors  Z1,/2, . . . ,/”  such 
that 

ZV  =  6ij  (9) 

for  i,j  =  1, . .  .  ,  n.  Here,  Z>,j  is  the  Ivronecker  delta  symbol: 

1  if  i  =  j, 

0  if  i  ^  j 

for  i,j  =  1, . . .  ,  n.  Eq.  (9)  will  be  referred  to  as  the  biorthonormality  property  of  the  vectors  Z* ,  rJ . 

Proposition  2  below  provides  a  convenient  representation  for  scalar  symmetric  multilinear  func¬ 
tions  on  ( IRn)k  in  terms  of  the  dual  basis  vectors  Z!  and  a  set  of  “structural  coefficients.”  The  result 
will  be  applied  in  the  next  section,  yielding  a  representation  of  Liapunov  function  candidates. 
Proposition  2.  (Coordinate  Representation  of  Scalar  Multilinear  Functions)  Any  sym¬ 
metric  fc-Iinear  function  ip  :  ( IRn)k  IR  can  be  written  as 

xp{x\x\  . . .  ,**)  =  V  ^.^(Z'VHZ’V) . . .  (Z****)  (11) 

tl  ,»2,...,i*  =  l 

where  the  (fc-tuple)  sum  is  taken  over  all  f i ,  *2 ,  ■  -  ■ ,  ,  and  where  the  structural  coefficients 

are  symmetric  with  respect  to  all  permutations  of  t'i,  *2,  •  •  • ,  ik- 

Proof.  This  result  is  a  special  case  of  [16,  Thm.  1.2].  However,  we  sketch  a  rather  straightforward 
proof  for  the  sake  of  completeness.  By  [4,  Thm.  2.12.2],  a  scalar  multilinear  function  is  determined 
by  its  values  when  evaluated  at  all  combinations  of  basis  vectors  as  arguments.  The  formula  (11) 
for  ip  is  clearly  that  of  a  (t-linear  function.  Moreover,  it  follows  from  the  biorthonormality  property 
(9)  of  basis  vectors  r*  and  dual  basis  vectors  Z’  that,  by  appropriate  assignment  of  the  structural 
coefficients  tpi1i2...ik,  any  set  of  such  values  may  be  achieved.  Hence,  the  representation  above  is 
sufficiently  general  to  accomodate  any  scalar  multilinear  function  ip.  □ 
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2.4  Complexification  of  Real  Multilinear  Functions 

It  is  sometimes  convenient  to  evaluate  a  (real)  multilinear  function  ip  :  ( IRn)k  i— ►  IRm  for  argument 
vectors  in  (Fn.  This  is  done  simply  by  evaluating  the  value  of  ip  as  if  the  argument  vectors  were  in 
IRn ,  using  a  representation  of  xp  such  as  Eq.  (11)  above.  This  process  is  the  complexification  of  ip. 

Our  use  of  the  complexification  device  is  relegated  to  Section  7,  in  the  construction  of  Liapunov 
functions  for  Eq.  (1)  under  hypothesis  (H).  The  following  observation  will  be  important  in  ensuring 
that  the  constructed  Liapunov  functions  are  indeed  real- valued. 

Proposition  3.  (Test  for  Realness  of  Multilinear  Functions)  Let  ip  denote  a  symmetric 
fc-linear  function  ip  :  ((T”)fc  i— »  (f771.  The  image  of  ( IRn)k  of  under  the  map  tp  is  IRm  if  and  only  if 

Hx\x2,...,xk)  -  ip{x\x2,...,xk)  (12) 

for  all  vectors  x1 ,  x2 , . . .  ,  xk  £  C" . 

Proof.  The  “if  ”  part  is  automatic.  An  induction  proof  is  now  sketched  for  the  “ only  if”  part. 
Let  j  denote  the  number  of  argument  vectors  x1  6  (Tn  that  are  not  also  in  IRn .  That  (12)  holds 
when  j  —  0  is  obvious.  Also,  if  (12)  holds  for  some  j  =  jo  <  k ,  then  it  is  a  simple  exercise  to  verify 
that  it  also  holds  for  j  —  jo  +  1.  □ 

3  REPRESENTATION  OF  LIAPUNOV  FUNCTION 
CANDIDATES 

Eq.  (1)  may  be  rewritten,  upon  Taylor  series  expansion  of  /(x),  in  the  form 

i  =  /(x) 

=  Lx  +  Q(x,  x)  +  C(x,  x,  x)  +  •  •  •  (13) 

Here,  L  :=  D/( 0)  and  Q(x,x),  C(x,x,x)  are  vector-valued  quadratic  and  cubic  forms,  with  the 
dots  denoting  higher  order  terms.  Without  loss  of  generality,  assume  that  Q(x,x)  is  induced  by 
a  symmetric  bilinear  function  Q(xl  ,x2)  and,  similarly,  that  C(x,x,x)  is  induced  by  a  symmetric 
trilinear  function  C(x1,x2,x3). 

In  Case  (S)  (one  zero  eigenvalue),  we  shall  in  the  sequel  seek  Liapunov  functions  V(x)  consisting 
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of  the  sum  of  a  quadratic  part  and  a  cubic  part,  viz. 

V(x)  =  xTVx  +  /C(x,  x,  x).  (14) 

It  is  natural  to  require  V  to  be  symmetric  and  positive  definite.  Similarly,  the  cubic  form  fC(x,  x,  x) 
is  induced  by  a  symmetric  trilinear  function  K{xx ,  x2,  x3).  Since  the  term  /C(x,x,x )  is  dominated 
by  the  quadratic  term  xTVx,  any  such  V(x)  will  indeed  be  locally  positive  definite. 

In  our  study  of  Case  (II)  (two  purely  imaginary  eigenvalues),  we  will  include  a  quartic  term 
T|x,x,x,x)  in  the  candidate  Liapunov  function  V(x)  in  addition  to  the  quadratic  and  cubic  terms 
presented  in  (2),  viz. 

V(ar)  =  xTVx  +  K(x,x,x) +  T(x,x,x,x).  (15) 

We  of  course  ask  that  the  quartic  form  T (x,  x,  x,  x)  be  induced  by  a  symmetric  tetralinear  function 
T(  x1 ,  x2,  x3,  x4).  Certainly,  the  local  positive  definiteness  of  V(x)  with  V  >  0  remains  preserved 
under  inclusion  of  T(x,  x,  x,  x)  or  terms  of  still  higher  order. 

Next,  we  invoke  Proposition  2  for  the  cases  k  =  2,3  and  4,  obtaining  coordinate  representations 
of  the  bilinear,  trilinear,  and  tetralinear  functions  xx^Vx2,  K{xx ,  x2,  x3),  and  T (x1 ,  x2,  x3,  x4), 
respectively. 

Consider  first  the  case  k  =  2.  Since  any  real  quadratic  form  xlTVx 2  is  determined  by  a  real 
symmetric  matrix  V,  we  can  apply  Proposition  2  to  conclude  that  all  such  matrices  V  have  the 
form 

V  =  (16) 

i=l  j=l 

where  n tj  =  irji  are  real  coeffiecients. 

The  following  representations  for  trilinear  functions  1C  and  tetralinear  functions  T  also  follow 
from  Proposition  2: 

/c(x\x2,x3)  =  (17) 

t=i  i=i  k= i 

T(x1,x2,x3,x4)  =  E  E  E  E  r>jfcp(/'x1)(Fx2)(^x3)(/Px4),  (18) 

i=l  j—  1  k=  1  p=  1 

respectively.  These  multilinear  functions  are  rendered  symmetric  by  imposing  the  condition  that 
the  values  of  the  structural  coefficients  (ntjk)  and  (rijkp)  do  not  depend  on  the  order  of  the  indices. 
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The  simple  representations  above  for  the  bilinear,  trilinear  and  tetralinear  functions  appearing  in 
the  Liapunov  function  candidates  V  imply  that  the  construction  ofV  is  tantamount  to  specification 
of  the  structural  coefficients  ir{j , 

In  the  calculations  to  follow,  the  goal  is  to  obtain  sets  of  structural  coefficients  7r ij,Kijk,rijkp 
which  ensure  the  local  negative  definiteness  of  V,  the  time  derivative  of  the  Liapunov  function 
candidate,  along  trajectories  of  (1).  Of  course,  this  will  only  be  possible  under  assumptions  on 
Eq.  (1)  which  guarantee  local  asymptotic  stability  of  the  origin.  In  the  next  section,  a  result  is 
presented  on  local  definiteness  of  a  class  of  bivariate  functions.  The  sufficient  conditions  for  local 
definiteness  provided  by  this  result  will  facilitate  a  systematic  derivation  of  local  Liapunov  functions 
and  conditions  for  local  asymptotic  stability. 

4  LOCAL  DEFINITENESS  OF  A  CLASS  OF 
BIVARIATE  FUNCTIONS 

We  now  introduce  an  interesting  lemma  which  will  prove  to  be  an  important  tool  in  exhibiting 
conditions  for  local  negative  definiteness  of  the  time  derivative  of  a  Liapunov  function  candidate. 
Lemma  1.  (Local  Definiteness  of  a  Class  of  Bivariate  Functions)  The  scalar  bivariate 
function 

6(u,v)  —  a2  on2  +  aoju4 

+  a2 \u2v  +  030 u3  +  «13 uu3  +  a22u2v2  +  o3iu3r  +  a40u4  +  0(|(u,  u)|5)  (19) 

in  the  real  variables  u  and  v  is  locally  negative  definite  near  (u,  v)  —  (0,  0)  provided  that  a2 o  <  0 
and  ao4  <  0.  Here,  0(|(u,u)|5)  denotes  terms  of  fifth  and  higher  order  in  |(u,u)|. 

Remark  1.  It  is  convenient  to  view  this  lemma  as  follows.  Consider  6(u,v),  a  scalar  polynomial 
function  of  the  scalars  u  and  v ,  for  which  the  leading  term  in  6(u,  0)  is  a2QU2 ,  and  the  leading  term 
in  6(0,  v)  is  ao4v4 .  Then  the  assertion  is  that  6(u,v)  is  locally  negative  definite  when  two  basic 
conditions  are  fulfilled:  First,  the  absence  of  the  two  terms  oq^v3 ,  a42uv2]  Second,  the  local  negative 
definiteness  of  the  univariate  functions  6(u,0)  and  6(0,  u)  (i.e. ,  a2o  <  0  and  ao4  <  0,  respectively). 
Proof.  In  proving  Lemma  1,  we  neglect  the  terms  0(|(u,u)|s)  in  6(u,v),  since,  being  higher  order 
terms,  they  can  easily  be  incorporated  with  only  slight  modifications  in  the  analysis.  With  this 
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understanding,  rewrite  6(u,v)  in  the  form  of  a  quadratic  polynomial  in  u: 


6(u,  v) 


(a2o  +  a21v  +  a30u  +  a22v 2  +  a3iuv  +  a40u2)u 2  +  a13v3u  +  a04v4 
p(u,  v)u2  +  (ai3u3)u  +  a04v4. 


(20) 


Here, 


p(u,  v )  :=  a 20  +  a24v  +  a30u  +  a22v2  +  a3im;  +  a40u2 . 

Since  p(0,  0)  =  a2o  <  0,  it  is  clear  that  there  is  an  >  0  such  that  p(u,  v)  <  0  for 


(21) 


M,  M  <  ei. 


(22) 


(One  could  easily  write  a  formula  for  such  an  fj.) 

The  leading  coefficient  p(u,v)  in  the  expression  (20)  for  6(u,v)  is  therefore  strictly  negative  for 
|u|,  |u|  <  ei .  Next,  rewrite  6(u,v)  as 


6(u,  v)  =  p(u,  d)[u  + 


ai3V3  -,2 
2  p(u,vy 


+  9(UT), 


(23) 


where 


4 p(u, v) 


{4020004 


+  4(130  a04u  +  4a2ia04v 
+  4a04a4ou2  ■+ 4a04a3iuv  +  (4a04a22 


«?3)^2}. 


(24) 


Since  a20  <  0  and  o40  <  0,  the  constant  term  in  the  expression  in  braces  in  Eq.  (24),  namely 
4020004,  is  strictly  positive.  Hence,  there  is  an  €2  >  0  such  that  the  expression  in  braces  in  Eq.  (24) 
is  strictly  positive  for 


M,M<e2.  (25) 

Recalling  that  p(u,  v)  <  0  for  |u|,M  <  tj,  we  have  that  for  |u|,  |n|  <  e  :=  min(fi,e2),  9(u,u)  <  0, 
with  q(u,  n)  =  0  only  for  v  =  0  (see  Eq.  (24)).  Now  consider  the  implications  of  these  observations 
for  the  expression  (23)  for  6(u,v).  Clearly,  for  |u|,  |u|  <  e  and  v  ^  0,  the  fact  that  q(u,v)  is  strictly 
negative  ensures  that  6(u,v)  <  0.  If,  on  the  other  hand,  |u|,  |v|  <  e  and  v  =  0,  then  6(u,v)  reduces 
to 


6(u,0)  =  p(u,0)u2 

<  0 


(26) 


for  u  0.  Thus,  6(u,v)  is  indeed  locally  negative  definite  near  (0,0). 


□ 
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5  CALCULATIONS  INVOLVING  THE  STABLE  SUBSPACE 


In  this  section,  we  define  the  stable  subspace  of  IRn  corresponding  to  the  Jacobian  matrix  L,  recall 
an  associated  orthogonality  property  from  [10],  and  employ  the  stable  subspace  concept  in  the 
choice  of  the  quadratic  term  xTVx  in  the  Liapunov  function  candidate  V(x)  (cf.  Eqs.  (14),  (15)). 
The  development  proceeds  for  Cases  (S)  and  (H)  in  parallel. 

Definition  4.  The  stable  subspace  of  IRn ,  denoted  by  Es ,  is  the  span  of  the  eigenvectors  (and 
generalized  eigenvectors,  if  any)  of  L  corresponding  to  the  stable  eigenvalues  of  L.  o 

In  Case  S,  any  vector  x  £  IRn  has  a  unique  representation  x  =  ar-{-w  where  a  is  a  real  scalar,  r 
is  the  right  eigenvector  of  L  corresponding  to  the  eigenvalue  0,  and  w  £  Es .  In  Case  H,  any  vector 
x  £  IRn  has  a  unique  representation  x  —  ar  +  ar  +  w  where  a  is  a  complex  scalar,  r  is  the  right 
eigenvector  of  L  corresponding  to  the  eigenvalue  iujc,  and  w  £  Es . 

The  following  property  is  well  known  (see,  e.g.,  [10,  Appendix  4.1]). 

Proposition  4.  (Orthogonality  of  Left  and  Right  Eigenvectors)  Let  P  and  denote  left 
and  right  eigenvectors,  respectively,  corresponding  to  eigenvalues  \a  and  of  a  matrix  A.  Either, 
or  both,  of  P  and  may  be  generalized  eigenvectors.  If  Xa  yf  \p,  then  Pr'9  =  0.  Moreover, 
the  subspace  of  all  column  vectors  nullified  by  P  is  precisely  the  span  of  all  right  eigenvectors  and 
generalized  right  eigenvectors  of  A  associated  with  eigenvalues  other  than  XQ. 

Remark  2.  Proposition  4  implies  the  following  facts,  which  will  prove  useful  in  the  sequel.  As 
above,  let  l  denote  a  left  eigenvector  corresponding  to  the  critical  eigenvalue  0  (in  Case  (S))  or  iu>c 
(in  Case  (H)).  Then  Iw  =  0  if  and  only  if  w  £  Es .  Moreover,  pw  =  0  for  a  row  vector  p  if  and 
only  if  p  £  span(/)  (Case  (S)),  or  p  £  span(5?e  /,9m  /)  (Case  (H)).  Finally,  in  Case  H,  we  have  Ir 
=  Jr  =  0. 

Since  the  Jacobian  matrix  L  =  Df{ 0)  has  part  of  its  spectrum  on  the  imaginary  axis,  it  is  not 
possible  to  choose  a  positive  definite  V  for  which  LT V  +VL  is  negative  definite.  However,  one  can 
ensure  that  the  latter  matrix  is  negative  definite  on  a  subspace  of  IRn,  while  being  only  negative 
semidefinite  on  all  of  lRn .  A  method  for  achieving  this  is  given  next. 

Recall  that  r  denotes  the  eigenvector  of  L  corresponding  to  the  critical  eigenvalue  (0  in  Case 
S,  iu>c  in  Case  (H)).  Note  that  r  £  <Vn  in  Case  H,  and  that  r  £  IRn  in  Case  (S).  Denote  by  Es  the 
subspace  of  IRn  spanned  by  the  eigenvectors  (and  generalized  eigenvectors,  if  any)  corresponding 
to  the  stable  eigenvalues  of  L  (in  either  Case  (S)  or  Case  (H)).  We  refer  to  Es  as  the  stable  subspace 
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of  lRn. 


The  following  proposition  is  useful  in  selecting  the  quadratic  term  xTVx  in  the  Liapunov  func¬ 
tion  V  under  either  hypothesis  (S)  or  (H). 

Proposition  5.  (Liapunov  Matrix  Equation  on  Stable  Subspace)  Using  the  notation  above, 
and  under  either  hypothesis  (S)  or  (II),  there  exists  a  family  of  real  symmetric  n  x  n  matrices  II 
for  which 

(i)  nr  =  0,  (ii)  wTUw  >  0,  and  (iii)  wT{LTU  +  UL)w  <  0  (27) 


for  all  w  £  Es ,  w  ^  0. 

Proof.  The  proof  may  be  carried  out  in  two  steps.  In  the  first  step,  we  exhibit  a  choice  of 
coordinates  for  the  state  space  IRn  for  which  the  existence  of  matrices  n  is  transparent,.  In  the 
second  step,  we  verify  that  the  existence  of  a  matrix  n  satisfying  (i)-(iii)  in  one  coordinate  system 
implies  the  existence  of  such  a  matrix  for  any  choice  of  coordinates.  Step  1.  Suppose,  then,  that 
the  state  x  of  system  (1)  is  expressed  with  respect  to  a  coordinate  basis  {r1 ,  r2, . . . ,  rn}  defined 
as  follows.  In  Case  S,  take  r1  r  and  choose  the  remaining  basis  vectors  r'  £  IRn,  i  —  2, . . .  ,n 
such  that  span{r2, .  . . ,  r"}  =  Es .  Analogously,  in  Case  H,  choose  r1  :=  Sfte(r),  r2  :=  3m(r),  and 
let  r3, . . . ,  r"  satisfy  span{r3, . . . ,  r”}  =  E* .  For  such  a  coordinate  basis,  L  has  the  block  diagonal 
representation 


L  - 


(28) 


Here,  Ls  is  a  real  square  stable  matrix  whose  eigenvalues  coincide  with  the  stable  eigenvalues  of  L, 


and  6  is  given  by 


f  0 

in  Case  S, 

( 0  Wc) 

in  Case  H. 

(29) 

K\-Uc  0  / 

(Note  that  9  is  a  scalar  for  Case  S,  and  is  a  2  x  2  matrix  for  Case  (H).)  It  is  now  straightforward 
to  exhibit  a  matrix  II  satisfying  (i)-(iii).  Consider  the  matrix 


(° 

\o  iw 


(30) 


where  n22  is  a  real  symmetric  positive  definite  matrix  of  dimension  (n—  1)  in  Case  S,  and  dimension 
(n  —  2)  in  Case  H,  for  which  the  matrix 


F  ,s^  II 22  +  n  22^5 


(31) 
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is  negative  definite.  The  existence  of  such  a  matrix  II22  is  clear,  since  La  is  stable.  Note  that,  for 
the  present  choice  of  coordinate  basis,  we  have  that  in  Case  (S)  the  first  component  of  any  vector 
w  £  Es  is  0,  and,  in  Case  H,  the  first  two  components  of  w  are  0.  Also,  the  right  eigenvector  r 
is  given,  in  Case  (S)  and  Case  H,  by  r  =  (1,  0, . . . ,  0)T  and  r  =  (1,  i,  0, . . . ,  0)r,  respectively.  The 
matrix  IT  of  Eq.  (30)  is  now  easily  verified  to  satisfy  conditions  (i)-(iii).  Step  2.  Next,  we  show 
that  existence  of  a  matrix  IT  satisfying  (i)-(iii)  in  one  set  of  coordinates  implies  existence  of  such 
a  matrix  for  any  set  of  coordinates.  Let  the  coordinate  change  be  determined  by  a  nonsingular 
transformation  matrix  E,  and  indicate  vectors  in  the  transformed  coordinates  by  a  hat  C): 

r  -  Er,  w  —  Hw.  (32) 

In  the  new  coordinates,  L  has  the  representation 

L  =  ELE-1.  (33) 

It  is  straightforward  to  verify  that  the  matrix 

II  =  (Er)-1nE“1  (34) 

satisfies  conditions  (i)-(iii),  where  all  quantities  are  taken  in  the  new  coordinates.  (Note  that  the 
relationship  between  II  and  II  specified  in  Eq.  (34)  is  not  that  of  a  similarity  transformation.)  o 

As  a  corollary  to  Proposition  5,  we  have  the  following  result. 

Corollary  1.  Let  IT  be  a  matrix  satisfying  (i)-(iii)  of  Proposition  5.  Then,  there  is  an  a  >  0  such 
that,  for  each  nonzero  w  £  Es, 

u>T(LT IT  +  II L)w  <  — a|u;|2.  (35) 

Proof.  This  follows  easily  using  the  proof  of  Proposition  5  and  the  standard  fact  that,  for  a 
negative  definite  real  matrix  Q ,  there  is  an  a  >  0  such  that  xTQx  <  — a|r|2  for  all  x  £  lRn .  In  the 
present  setting,  the  role  of  Q  is  played  by  LST II22  +  IT 22 Ls-  D 

In  our  construction  of  Liapunov  functions  for  Case  S,  we  shall  employ  matrices  V  of  the  form 

v  =  n  +  iTi  (3C) 
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where  II  is  any  real  nxn  matrix  satisfying  the  conditions  of  Proposition  5.  For  Case  II,  the  matrices 
V  will  be  of  the  form 

V  =  U  +  lTl+JTl.  (37) 

It  is  not  difficult  to  check  that,  both  in  Case  (S)  and  Case  H,  the  matrix  V  is  positive  definite,  and 
the  matrix  L? V  -\-VL  is  negative  semidefinite. 

6  CONSTRUCTION  OF  LIAPUNOV  FUNCTIONS  IN  THE 
CASE  OF  ONE  ZERO  EIGENVALUE 

In  this  section,  we  construct  a  family  of  Liapunov  function  candidates  for  system  (1)  (equivalently, 
(13))  under  hypothesis  (S).  The  main  task  is  to  specify  the  matrix  V  and  the  cubic  form  K( x,x,x) 
appearing  in  the  expression  (3)  for  V(x).  In  the  foregoing,  we  have  constrained  the  matrix  V  to 
take  the  form  V  —  lTl  +  II  in  Case  S,  where  II  is  any  real  symmetric  nxn  matrix  satisfying  the 
conditions  of  Proposition  5. 

6.1  Conditions  for  V  <  0 

Using  Proposition  1  (Jacobian  of  Homogeneous  Functions),  the  time  derivative  of  V(x)  evaluated 
along  trajectories  of  Eq.  (13)  can  be  written  as 

V(x)  =  xT(LTV  +  VL)x 

+  '2Qt(x,x)Vx  +  ZK(x,  x,  Lx) 

+  2 CT(x,x,  x)Vx  +  3 /C(x,  x,  Q(x,  x))  +  •  •  •  (38) 

Recall  from  Section  5  the  notation  Es  for  the  stable  subspace  of  IRn ,  i.e.,  the  (n—  l)-dimensional 
subspace  spanned  by  the  eigenvectors  (and  generalized  eigenvectors,  if  any)  corresponding  to  the 
stable  eigenvalues  of  L.  Using  the  representation  x  —  ar  +  ui  ( w  E  Es)  in  Eq.  (38);  recalling  that 
P  has  been  chosen,  by  Eq.  (36),  such  that  V  =  ITl  -f  n  with  nr  =  0  ;  invoking  the  fact  that  a 
multilinear  function  is  linear  in  each  argument;  and  collecting  terms  on  the  right  side  of  (38)  of  like 
order  in  |(a,  ™)l>  we  obtain  a  series  expansion 

v(x)  [y](2)  +  [y](3j  +  [y](4)  +  •  •  •  (39) 
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where  the  integer  subscripts  denote  the  degree  of  the  corresponding  term  in  |(a,in)|,  and  the  dots 
denote  terms  of  fifth  and  higher  order  in  |(a,  tw)|.  Specifically,  the  terms  appearing  on  the  right  side 
of  Eq.  (39)  are  given  by 

[V](2)  =  wT(LTn  +  UL)w,  (40) 

P>](3)  =  2a3IQ(r,  r) 

+  3a2Af(r,  r,  Lw)  4-  4 a2lQ{r,  w )  +  2 a2QT(r,  r)Ilu> 

+  2alQ{u>,  w)  +  4aQT(r,  in)IIu>  +  6a£(r,  w,  Lw) 

+  2  QT  (w,w)Nw  4- 3E(w,w,  Lw),  (41) 

M(4)  =  a4{2/C(r,r,  r)  +  3AJ(r,r,  Q(r,  r))} 

4-  2a3{3/C(r,  r,  w)  4-  3 K(r,  w,  Q(r,  r))  +  3 fC(r,  r,  Q(r,  in))} 

4-  3a2{2/C(r,  w,  w)  +  /C(r,  r,  Q(w,  w))  +  4 K{r,  w,  Q(r,  w))  +  K(w,  w,  Q(r,  r))} 

4-  2a{lC(w,  w,  w)  4-  3 fC(r,  w,  Q(w,  in))  4-  3/C(w,  w ,  Q(r,  in))} 

4-  2CT(w,  in,  in)IIu>  4-  3AC(tn,  in,  Q(w,  in)).  (42) 

Note  that  condition  (iii)  on  II  (cf.  Proposition  5)  implies  that  [V]^2j  <  0  for  in  £  Es,w  ^  0. 

This  does  not  of  course  imply  that  V  is  locally  negative  definite,  only  that  it  is  locally  negative 

definite  on  the  subspace  Es . 

Lemma  1,  along  with  the  foregoing  computation  of  V,  allow  us  to  obtain  the  following  pre¬ 
liminary  statement  concerning  the  local  asymptotic  stability  of  the  origin  of  Eq.  (1).  Note  that 
condition  (SI)  of  the  next  proposition  is  a  known  necessary  condition  for  stability  for  systems  (1) 
possessing  a  simple  zero  eigenvalue  (see  for  instance  [2],  [10]). 

Proposition  6.  Under  hypothesis  (S),  the  origin  of  Eq.  (1)  is  locally  asymptotically  stable  if 
there  are  a  real  symmetric  nx  n  matrix  II  satisfying  (i)-(iii)  of  Proposition  5,  and  a  symmetric  real 
trilinear  function  A^x1 ,  x2,  x3),  for  which  the  following  three  conditions  hold: 

(51)  lQ(r,r)  =  0, 

(52)  3 fC(r,  r,  Lw)  4-  4 lQ(r,  w)  4-  2 QT(r,  r)IIu;  =  0  for  all  w  £  Es ,  and 

(53)  2/C(r,  r,  r)  4-  3A f(r,  r,  Q(r,  r))  <  0. 
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Proof.  Let  conditions  (S1)-(S3)  of  the  Proposition  hold.  It  is  straightforward  to  write  an  upper 
bound  for  V(ar  +  w )  in  the  form  of  a  scalar  bivariate  Taylor  series.  To  facilitate  application  of 
Lemma  1,  we  employ  notation  consistent  with  Eq.  (19)  of  Lemma  1,  and  define  variables  u  :=  |u;| 
and  v  :=  |a|.  The  proof  proceeds  in  two  steps.  The  first  step  consists  of  verifying  that  093  =  0  and 
ai2  =  0.  In  the  second  step,  we  ascertain  that  a 20  <  0  and  ao4  <  0.  Step  1.  From  Eq.  (41),  it  is 
clear  that  the  u3-term  in  V  is  2 lQ(r,r)v3,  and  this  vanishes  by  virtue  of  condition  (SI).  Thus,  the 
upper  bound  for  V  will  naturally  be  absent  of  a  term  a 03U3.  Using  Eq.  (41),  it  is  apparent  that 
the  uu2-term  in  the  upper  bound  will  vanish  if  3 A(r,  r,  Lw)  +  4/Q(r,  w)  +  2 Q1  (r,  r)IIu;  =  0  for  each 
w  G  Es.  This  latter  condition  is  precisely  (S2).  Step  2.  Eq.  (40)  and  Coro. ..'ary  1  imply  that  the 
quadratic  terms  in  V  are  bounded  above  by  a  function  a2ou2,  where  a2o  <  0.  Also,  Eq.  (42)  (which 
gives  the  quartic  terms  in  V)  and  assumption  (S3)  together  imply  that  a04  <0.  d 

6.2  Algorithm  for  Construction  of  V  in  Case  (S) 

For  Proposition  6  to  be  useful  in  the  explicit  construction  of  Liapunov  functions  for  Eq.  (1),  a 
method  is  needed  facilitating  the  choice  of  a  matrix  II  and  a  trilinear  function  K  for  which,  under 
an  auxiliary  condition  guaranteeing  local  stability,  (S2)  and  (S3)  are  satisfied.  As  it  turns  out,  one 
can  first  choose  any  matrix  II  for  which  conditions  (i)-(iii)  of  Proposition  5  hold,  and  then  proceed 
to  construct  compatible  trilinear  functions  K  satisfying  (S2),  (S3).  We  now  proceed  to  construct  a 
family  of  such  trilinear  functions,  using  the  representation  (17)  in  terms  of  the  associated  structural 
coefficients. 

The  general  representation  (17)  for  trilinear  functions  Af(x] ,  x2 ,  x3)  assumes  a  specific  choice  of 
basis  for  IRn .  Let  {r1 ,  r2,  . . .  ,  rn }  be  a  basis  for  IRn ,  obtained  by  setting  r1  :=  r  and  requiring 
that  rl  £  Es  for  i  =  2 Let  the  associated  dual  basis  (discussed  in  Section  2.3)  be  given 
by  {l1 ,  /2, . . . ,  ln}.  Recall  that  the  biorthonormality  property  (9)  holds,  i.e.,  that  Pr3  =  6tj ,  the 
Kronecker  delta.  Moreover,  by  Proposition  4  and  the  fact  that  the  dual  basis  is  unique,  we  have 
that  f1  =  /,  the  left  eigenvector  of  L  associated  with  the  eigenvalue  zero  (recall  the  normalization 

fr=l)- 

By  (17),  the  trilinear  function  K  can  be  represented  as  follows: 

K(x\x2,x3)  =  (43) 

t=l  j=l  k=  1 
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By  (SI)  and  Proposition  4,  we  have  that  Q(r,r )  €  Es .  Thus  we  are  free  to  substitute  Q(r,r ) 
for  Lw  in  condition  (S2),  upon  which  we  directly  obtain 

3 lC(r,r,Q(r,r))  -  -AlQ{r,  L~ Q{r,r))  -  2QT (r,r)UL~ Q(r,r) 

=  -4lQ(r,L-Q(r,r))-QT(r,r)((L~)TU +  UL-)Q(r,r).  (48) 

Here,  L~  is  as  defined  in  Eq.  (45).  Substitution  of  (48)  into  (S3)  yields  the  following  condition 
equivalent  to  (S3): 

2 lC(r,  r,  r)  -  4/Q(r,  L~Q(r,  r))  -  QT(r,  r)((L~f  n  +  H L~)Q(r,  r )  <  0.  (49) 

Note  that  (S3)  therefore  places  a  condition  only  on  the  quadratic  part  of  the  Liapunov  function 
candidate  V,  as  reflected  by  the  appearance  in  (49)  of  the  matrix  II. 

Conditions  (SI)  and  (S3)  are  akin  to  conditions  that  arise  in  the  stability  analysis  of  stationary 
bifurcation  for  parametrized  embeddings  of  Eq.  (1).  Under  hypothesis  (S),  such  parametrized 
systems  will  generically  exhibit  a  bifurcation  in  which  a  new  equilibrium  xt  coexists  with  the  origin 
for  each  small  |e|.  Here,  e  is  a  (normalized)  real  amplitude  parameter.  The  eigenvalue  near  zero  of 
the  bifurcated  equilibrium  xc  is  given  by  an  expansion 


/?(e)  =  /?1e  +  /?2e2  +  /?3e3  +  ---  (50) 

To  guarantee  asymptotic  stability  of  the  new  equilibrium  x£,  one  requires  —  0  and  /?2  <  0.  In 
this  context  (cf.  [2]),  we  have  the  “bifurcation  formulae” 


0\  =  lQ{r,r), 


(51) 


h  =  2/{C(r,  r,  r)  -  2 Q(r,  L~Q(r,  r))}.  (52) 

Note  that  condition  (SI)  is  therefore  identical  to  the  bifurcation  stability  condition  f3\  =  0. 
Similarly,  (S3)  (equivalently,  Eq.  (49))  is  readily  expressed  in  terms  of  the  coefficient  02-  Denote 
by  As(n)  the  n-dependent  scalar 

As(II)  :=  - QT(r ,  r)((Zr)Tn  +  nT~)Q(r,  r).  (53) 
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Eq.  (49)  is  now  rewritten  in  terms  of  fa  and  A5: 

fa  +  As(II)  <  0.  (54) 

The  following  result  is  well  known  (see,  e.g.,  [2],  [9]). 

Theorem  1.  Let  hypothesis  (S)  hold,  and  suppose  that  fa  =  0  and  fa  <  0,  where  fa  and  fa  are 
given  by  (51)  and  (52),  respectively.  Then  the  origin  of  Eq.  (1)  is  locally  asymptotically  stable. 

In  our  pursuit  of  Liapunov  functions  for  (1),  we  have  in  fact  rederived  this  result.  Indeed,  note 
that  our  condition  (SI)  requires  that  fa  =  0,  and,  for  a  given  system  (1)  for  which  ffa  <  0,  the 
matrix  II  can  be  chosen  so  as  to  ensure  that  (54)  holds.  (Recall  that  II  is  any  real  symmetric  matrix 
satisfying  conditions  (i)-(iii)  of  Proposition  5.  These  conditions  are  linear  in  II.) 

An  observation  of  relevance  here  is  that,  for  any  choice  of  II,  the  quantity  As(II)  is  nonnegative: 

As(n)  >  0.  (55) 

Thus,  As(II)  contributes  adversely  to  satisfaction  of  Eq.  (53). 

We  can  now  present  the  main  result  of  this  section,  the  construction  of  a  family  of  Liapunov 
functions  V(z)  of  the  form  (14)  for  Case  (S)  (one  zero  eigenvalue). 

Theorem  2.  Let  hypothesis  (S)  hold,  and  suppose  that  fa  =  0  and  fa  <  0,  where  fa  and  fa  are 
given  by  (51)  and  (52),  respectively.  Then  any  function  V(x)  resulting  from  Algorithm  Vs  below 
is  a  Liapunov  function  for  the  equilibrium  point  0  of  Eq.  (1). 

Algorithm  V5.  (Construction  of  Liapunov  functions  V(x)  =  xTVx  +  IC(x,x,  x )  in  Case  (S)) 

Step  1.  Compute  l  and  r.  Choose  a  basis  {r2,...,r"}  for  Es .  Compute  the  dual  basis 

{l1,  /2, . .  . ,  /"}  to  the  basis  {r1 ,  r2, . . . ,  r”}  for  lRn .  Here,  r1  :=  r  and  / 1  =  /.  Compute 

the  coefficients  fa  and  fa  according  to  Eqs.  (51)  and  (52),  respectively.  Check  that  fa  =  0 
and  fa  <  0. 

Step  2.  Choose  any  real  symmetric  n  x  n  matrix  n  satisfying,  for  all  w  E  Es ,  w  7^  0:  (i)  Hr  =  0, 

(ii)  ujtUw  >  0,  and  (iii)  wT(LTU  +  IU)®  <  0,  and  for  which 

|As(n)|  <  \fa\,  (56) 

where  As(n)  is  as  defined  in  Eq.  (53). 
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Step  3.  Set  V  =  lTl+  IL 


Step  4.  Set  the  structural  coefficients  K;n,  i  =  2, . . . ,  n  to 

K«n  =  ”{2 lQ(r,L~r')  +  QT(r,  r)UL~r'}.  (57) 

Step  5.  Symmetry  requires  that  be  independent  of  permutations  in  the  indices  The 

structural  coefficients  in  the  representation 

IC(x\x2,x3)  =  EEE  Kut(/’xl)(/^2)(/fcx3)  (58) 

t— 1  j— 1  k=  1 

which  have  not  been  specfied  in  Steps  1-4  are  either  determined  by  symmetry  and  Eq.  (57), 
or  can  be  assigned  arbitrarily,  subject  only  to  the  symmetry  requirement. 

7  CONSTRUCTION  OF  LIAPUNOV  FUNCTIONS  IN  THE 
CASE  OF  A  PAIR  OF  PURE  IMAGINARY  EIGENVALUES 

In  this  section,  Liapunov  functions  are  constructed  for  the  origin  of  Eq.  (1)  under  hypothesis  (II). 
The  construction  parallels  that  of  the  last  section,  while  differing  from  it  in  several  respects.  For 
example,  the  Liapunov  function  candidates  used  in  the  previous  section  consist  only  of  quadratic 
and  cubic  terms  in  the  state.  In  this  section,  quartic  terms  also  appear  in  the  Liapunov  function 
candidates  (cf.  Eq.  (15)).  If  quartic  terms  were  not  included  in  the  assumed  form  of  the  Liapunov 
function  candidates,  the  construction  would  fail  generically  in  Case  (H).  A  second  difference  be¬ 
tween  the  calculations  of  this  and  the  preceding  section  concerns  the  adoption  of  complex  notation. 
Although  not  essential,  complex  notation  is  both  natural  and  convenient  when  considering  local 
asymptotic  stability  of  systems  (1)  under  hypothesis  (H). 
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7.1  Conditions  for  V  <  0 


Consider,  then,  the  time-derivative  of  the  Liapunov  function  candidate 

V(x)  =  xTVx  +  K(x,  x,  x)  -f  T(x,  x,  x,  x)  (59) 

along  trajectories  of  Eq.  (1).  Using  Proposition  1,  and  the  Taylor  series  representation  (13)  of 
/(x),  we  find  that  this  derivative  is  given  by 

V(x)  =  xT{LTV  +  VL)x 

+  “2QT(x,  x)Vx  +  3£(x,  x,  Lx) 

+  2 CT(x,  x,  x)Vx  +  3 K(x,  x,  Q(x,  x))  +  4T (x,  x,  x,  Lx)  +  •  •  ■  (GO) 

In  Section  5,  it  was  noted  that  any  vector  x  £  IRn  has  a  unique  representation 

x  =  ar  +  ar+w.  (61) 

Here,  r  is  the  right  eigenvector  of  L  associated  with  the  eigenvalue  iuic  (as  specified  in  Section  1),  r 
is  the  complex  conjugate  of  r  and  as  such  is  a  right  eigenvector  of  L  associated  with  the  eigenvalue 
—  iurc,  a  is  a  complex  scalar,  and  w  £  Es . 

The  first  step  in  our  procedure  for  obtaining  conditions  for  local  negative  definiteness  of  V(x) 
is  to  substitute  in  Eq.  (60)  the  representation  (61)  for  x  and  the  representation  (37)  for  the  matrix 
V,  and  group  terms  according  to  their  order  in  |(a,d,m)|.  (Recall  that  (37)  states  that,  in  Case 
II,  V  is  chosen  from  among  matrices  of  the  form  V  =  n  +  l1 1  +  /  TI,  where  n  satisfies  conditions 
(i)-(iii)  of  Proposition  5.)  Using  Proposition  4  and  Remark  2  (orthogonality  of  left  and  right 
eigenvectors),  and  after  a  considerable  amount  of  algebra  and  reordering  of  terms,  we  can  obtain 
explicit  formulae  for  the  quadratic,  cubic  and  quartic  terms  in  |(a,a,u;)|  (these  are  [V](2p  Mw 
and  [V](4j,  respectively)  in  the  expansion 

V(x)  :=  [V](2)  +  M(3)  +  1X1(4)  7 -  (62) 

of  V(x).  (This  is  in  analogy  with  Eqs.  (39)-(42)  of  the  preceding  section.)  For  example,  one  can 
check  that  [V](2p  Mp)  are  given  by 

[V](2)  =  wT(LTn  +  UL)w,  (63) 
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[V](3)  =  a3{2lQ(r,r)  +  3iwcJC(r,r,r)} 

+  a3{2lQ(r,  r)  -  3iuic/C(r,  r,  f)} 

+  a2a{2lQ(r,  r )  +  4/Q(r,  r)  +  3iwcfC(r,  r,  r)} 

+  a2a{2lQ(r,  r)  +  4/Q(r,  r)  -  3 iu>cfC(r,  r,  r)} 

+  a2{2QT(r,  r)IItD  +  4 /Q(r,  id)  +  3C(r,  r,  Lid)  +  6iu>cK(r,  r,  id)} 

+  a2{2QT(r,  r)IIui  +  4 lQ(r,  w)  +  3 IC(f,  r,  Lw )  —  6 iu>cK-(r,  f,  id)} 

+  aa{4QT(r,  r)IIiD  +  4/Q(f,  w)  +  4/Q(r,  id)  +  6AT(r,  f ,  Ltd)} 

+  a{2/Q(tD,  id)  +  QfC(r,  w,  Lw)  +  2>iwcK{r,  w,  id)} 

+  a{2lQ(w,  w)  -f  6£(r,  w,  Lw)  —  2>iwcK(f,w,  id)} 

+  3IC(w,w,  Lw).  (64) 

Note  that  both  [V]j2)  and  [V](3j  as  given  above  are  real-valued,  as  expected.  There  is  no  need  to 
give  the  full  expression  for  [V]^4j  here.  Instead,  we  proceed  directly  to  the  statement  of  the  following 
preliminary  result  in  the  construction  of  Liapunov  functions  under  hypothesis  (H).  In  the  proof  of 
this  result,  which  is  analogous  to  Proposition  6  in  the  preceding  section,  values  of  certain  pertinent 
terms  appearing  in  the  expansion  of  [V](4)  will  be  given.  All  the  terms  in  the  expansion  of  [V]^4j 
may  be  obtained  readily,  to  result,  in  an  expression  analogous  to  (64). 

Proposition  7.  Let  hypothesis  (H)  hold.  Suppose  that,  for  some  real  matrix  II  satisfying  the 
conditions  of  Proposition  5,  symmetric  real  trilinear  function  K\xl ,  x2 ,  x3),  and  symmetric  real 
tetralinear  function  T(x1 ,  x2,  x3,  x4),  the  following  seven  conditions  hold: 

(HI)  2 lQ(r,  r)  +  3 iwc/C(r,  r,  r)  =  0, 

(H2)  2 lQ(r,  r)  +  4 JQ(r,  r )  +  3iwc£(r,  r,  r)  =  0, 

(H3)  2 Q1  (r,  r)IIin  +  4lQ(r,  w)  +  3 fC(r,  r,  (L  +  2 iu>cI)w)  =  0,  for  all  w  E  Es , 

(H4)  2Q  J  (r,r)Ilw  +  2lQ(r,w)  +  2lQ(r,w)  +  3tC(r,r,  Lw)  =  0,  for  all  w  E  Es , 

(H5)  2 lC{r,  r,  r)  +  3 K{r,  r,  Q(r,  r))  -f  4iwcT (r,  r,  r,  r)  =  0, 

(HG)  lC(r,  r,  r)  +  3 fC(r,  r,  r)  +  +3A C(r,  r,  Q(r,  r))  +  3 K(r,  r,  Q(r,  r))  +  4iwcT (r,  r,  r,  r)  =  0, 
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(H7)  «Re{2 lC(r,  r,  r)  +  £(r,  r,  Q(r,  r))}  +  2 K(r,  r,  Q(r,  r))  <  0. 

Then  the  origin  of  Eq.  (1)  is  locally  asymptotically  stable. 

Proof.  The  proof  consists  of  a  judicious  application  of  Lemma  1,  which  gives  general  sufficient 
conditions  for  local  negative  definiteness  of  a  class  of  bivariate  functions.  Identify  the  scalar  variables 
u  and  v  in  Lemma  1  as  u  :=  |uj,  v  :=  |a|  =  |a|.  The  bivariate  function  <5(u,  v)  of  Lemma  1  is  taken 
to  be  a  local  upper  bound  on  V(z)-  We  show  that  conditions  (H1)-(II7)  are  sufficient  for  there  to 
exist  such  an  upper  bound  6(u,  v)  satisfying  the  hypotheses  of  Lemma  1.  To  ensure  absence  of  the 
u3-term  in  S(u,  v),  i.e. ,  that  a 03  =  0,  we  require  the  coefficients  of  a3,  a3,  a2a  and  d2a  in  Eq.  (64) 
to  vanish.  The  coefficient  of  a3  is  precisely  the  expression  on  the  left  side  of  (HI).  (Note  that  the 
a3-term  in  (64)  is  the  complex  conjugate  of  the  a3-term,  and  thus  also  vanishes  when  (III)  is  in 
force;  analogous  comments  apply  below.)  Similarly,  the  expression  on  the  left  side  of  (H2)  is  simply 
the  coefficient  of  a2a  in  Eq.  (64).  Thus,  (HI)  and  (H2)  combined  ensure  that  a 03  =  0.  To  ensure 
absence  of  the  uu2-term  in  6(u,v ),  i.e.,  that  a\2  =  0,  we  require  the  “linear-in-u>”  coefficients  of 
a 2  and  ad  in  Eq.  (64)  to  vanish  for  each  u>  g  Es.  (The  coefficient  of  a2,  being  the  conjugate  of 
that  of  a2,  will  then  vanish  automatically.)  Inspection  of  Eq.  (64)  reveals  that  this  is  equivalent  to 
conditions  (H3)  and  (H4)  above.  It  remains  to  show  that  <220  <  0  and  a40  <  0  (in  the  notation  of 
Lemma  1).  That  020  <  0  follows  immediately  from  Eq.  (63)  and  Corollary  1.  Conditions  ensuring 
that  a4o  <  0  can  only  result  from  examination  of  the  quartic  terms  in  V.  However,  quartic  terms  in 
[V](4)  which  involve  w  are  irrelevant  to  this  requirement.  The  needed  coefficients  can  be  obtained 
readily  by  substituting  Eq.  (61)  for  x  in  the  formula 

[V](4)(i)  =  2  CT(x,x,x)Vx  +  3  fC(x,x,Q(x,x))  +  4  T(x,x,x,Lx),  (65) 

using  the  orthogonality  of  left  and  right  eigenvectors  (cf.  Proposition  4  and  Remark  2),  Eind  using 
the  fact  that  a  multilinear  form  is  linear  in  each  argument.  The  expansion  of  [V](4)(ar  +  ar  -f  w) 
is  seen  to  contain  five  terms  which  are  quartic  in  (a,  a):  an  a4-term,  an  a3a-term,  the  conjugates 
of  these  two,  and  an  a2a2-term.  Of  these  five,  only  the  latter  is  sign-definite:  a2a2  =  |a|4.  Thus, 
we  require  the  coefficient  of  a2 a2  in  the  expansion  of  [V]^4)(ar  -f  ar  +  w)  to  be  negative,  and  the 
coefficients  of  a4  and  a3a  to  vanish.  It  is  readily  verified  that  the  left  side  of  (115)  is  the  coefficient 
of  a4,  and  the  left  side  of  (H6)  is  half  the  coefficient  of  a3d.  Finally,  it  is  straightforward  to  check 
that  the  left  side  of  (BI7)  is  one-sixth  the  coefficient  of  a2 a2  in  the  expansion  of  [V](4)(ar  +  ar  +  w). 
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Thus,  conditions  (H5)-(H7)  together  give  the  desired  negativity  of  ao4-  o 

7.2  Algorithm  for  Construction  of  V  in  Case  (H) 

Conditions  (H1)-(H7)  may  be  solved  for  a  trilinear  function  K(xl ,  x2,  x3)  and  a  tetralinear  function 
7  (x1  ,  x2,  x3,  x4),  under  an  appropriate  auxiliary  condition  ensuring  stability  of  the  origin  and  for 
a  given  matrix  II  satisfying  the  conditions  of  Proposition  5.  The  procedure  is  much  the  same  as 
was  carried  out  in  the  preceding  section,  where  the  coordinate  representation  of  a  trilinear  form 
K (x1 ,  x2 ,  x3)  was  employed  to  solve  (S1)-(S3)  for  the  structural  coefficients  Kijk  of  fC.  Due  to  this 
similarity,  only  a  summary  of  the  main  steps  in  the  derivation  is  deemed  necessary  here,  with  the 
result  for  the  Liapunov  functions  we  obtain  summarized  below  in  Algorithm  Vh ■ 

For  convenience,  we  continue  to  employ  complex  notation,  and  choose  a  coordinate  basis 
{  r1 ,  r2 , . .  . ,  r”  }  for  IRn  in  which  r1  :=  r,  r2  :=  r,  and  r3, . . . ,  rn  lie  in  Es  C  IRn  ■  To  this  ba¬ 
sis  there  corresponds  a  unique  dual  basis  of  row  vectors  {P ,  l2, . . . ,  /”},  where  P  :=  /  and  l2  :=  J. 
The  trilinear  and  tetralinear  functions  £(x] ,  x2 ,  x3)  and  T(x1 ,  x2,  x3,  x4)  are  then  expressed  in  the 
coordinate  representations  (17)  and  (18),  respectively.  We  seek  the  minimum  set  of  specifications  on 
the  associated  structural  coefficients  Ktjk  and  r,jkp,  respectively,  under  which  conditions  (II 1  )-(H 7) 
above  hold. 

Conditions  (HI)  and  (H2)  are  interpreted  in  this  framework  simply  as  assigning  values  to  the 
structural  coefficients  «ni  and  ^112>  respectively.  Next  consider  (H3).  Since  Es  is  invariant  under 
L,  w  £  Es  implies  that  (L  +  2iu>cI)w  lies  in  the  complexification  of  Es  (also  referred  to  as  Es  below). 
Thus,  we  can  define  a  vector  w  (L  +  2iu>cI)w,  noting  that  the  matrix  inverse  in  the  equation 
w  =  (Z,  +  <2wcl)~1w  exists  by  hypothesis  (H).  Interpreting  (H3)  as  a  requirement  on  each  basis 
vector  r3, . . . ,  r"  of  Es ,  we  find  that  (113)  amounts  to  a  specification  of  the  structural  coefficients 
Kill,  i  —  3,4,...,n.  Similarly,  (H4)  amounts  to  a  specification  of  the  structural  coefficients  km, 
i  —  3,4,...,  n. 

Since  each  of  the  structural  coefficients  Km,  i  =  l,2,...,n  is  fixed  (by  one  of  (H1)-(H3)), 
K{r,  r,  x)  is  fixed  for  any  x  £  (Pn .  The  coefficients  Km,  i  —  1, 2, . . . ,  n  are  also  fixed:  K112  is  fixed  by 
(H2),  Kii2,  i  =  3,4, . . .  ,n  are  fixed  by  (H4),  and  we  also  have  K212  =  KU2  by  applying  Proposition 
3  to  (112).  Thus,  /C(r,r,x)  is  also  fixed  for  any  x  £  6Jn .  By  these  remarks,  it  follows  that  the 
terms  AC(r,  r,  Q(r,  r)),  fC(r,  r,  Q(r,  f)),  K-(r,  f,  Q(r,  r))  appearing  in  (H5)  and  (H6)  are  determined  by 
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(H1)-(H4).  Their  values  may  be  found  by  expressing  Q(r,r )  and  Q(r,  r)  as  linear  combinations  of 
the  basis  vectors  r’,  i  =  and  then  employing  the  coordinate  representation  of  fC.  Thus, 

(H5)  and  (H6)  serve  to  assign  the  values  of  T(r,r,r,r)  (=  run)  and  T(r,  r,  r,  r )  (=  rm 2)-  The 
importance  of  including  the  quartic  term  T  in  the  Liapunov  function  candidate  now  becomes  clear: 
with  T  =  0,  (H5)  and  (H6)  become  constraints  on  the  system  which  do  not  constitute  necessary 
conditions  for  stability.  However,  with  inclusion  of  a  quartic  term  T,  (H5)  and  (H6)  are  quite  easily 
satisfied. 

By  the  remarks  above,  it  follows  that  the  quantity  appearing  on  the  left  side  of  (H7)  is  completely 
specified  by  the  system  dynamics  and  the  matrix  n.  Next,  we  sketch  the  derivation  of  an  explicit 
reformulation  of  the  left  side  of  (H7)  in  terms  of  system  (1)  and  n.  A  stability  coefficient  which 
arises  in  the  study  of  Hopf  bifurcation  for  parametrized  versions  of  (1)  under  hypothesis  (PI)  will 
appear  in  the  reformulation.  The  value  of  this  coefficient,  which  we  denote  as  /?2,  is  recalled  next. 
Note  that  this  coefficient  /32  is  distinct  from  the  coefficient  of  the  same  name  appearing  in  Section 
6.  In  the  present  context,  /?2  relates  to  an  (even  in  e)  expansion  /3(e)  =  /?2e2  +  /^c4  +  ■  •  •  of  the 
Floquet  exponent  near  zero  of  bifurcated  periodic  solutions  of  parametrized  embeddings  of  Eq.  (1). 
(Compare  with  Eq.  (50)  for  the  analogous  eigenvalue  expansion  in  Case  (S).) 

Define  vectors  £  and  77  by 

£  :=  -^-1Q(Dr),  (66) 

77  :=  ^(2iucI  -  L)_1Q(r,  r).  (67) 

(In  [9]  and  [1],  £  and  77  are  denoted  as  a  and  b ,  respectively.)  Then  /?2  is  given  by  the  “bifurcation 
formula”  [9],  [1] 

h  :=  23?e{2/Q(r,0  +  lQ(f,  77)  +  jj/C(r,  r,  r)}.  (68) 

Conditions  (H1)-(H4)  can  be  used  to  replace  (H7)  with  an  equivalent  condition  stated  explicitly 
in  terms  of  the  Taylor  expansion  (13)  of  system  (1).  Consider  the  tw’o  terms  fC(r,r,Q(f,f))  and 
K‘(r,  r,  Q(r,  r))  appearing  in  (H7).  The  former  quantity  is  of  the  form  fC(r ,  r,  x),  which  occurs  in  the 
statements  of  conditions  (HI)  (with  x  =  r ),  (H2)  (with  x  =  r)  and  (H3)  (with  x  €  Es).  Thus,  we 
resolve  the  vector  x  (=  Q(r,  r))  into  its  components  lying  in  the  subspace  Es  and  along  the  r -  and 
r-directions,  and  then  apply  (II1)-(H3).  The  r-component  of  any  vector  x  is  given  by  (/x)r  (recall 
the  normalization  Ir  =  1),  the  f-component  is  ( lx)r ,  and  the  ^-component  is  the  remainder 
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(69) 


xs  :=  x  —  (Ix'jr  —  (/x)r. 

Using  (H1)-(II3),  we  have  that,  for  any  x, 

K{r,r,x)  =  -JL-(lx)lQ(r,r)-  JL-(lx)(lQ(r,r) +  2lQ(r,f)) 

-  ^{2QT(r,  r)II(L  +  2iwcI)~lxs  +  4 lQ(r,  (L  +  2 «*/)"  V}.  (70) 

Since  r  is  an  eigenvector  of  L,  it  is  also  an  eigenvector  of  (L  +  2 iwcJ)_1  and  of  (L  -  2wc/)_1. 
Using  this  fact,  and  the  fact  that  fir  =  0,  this  expression  may  be  expanded  and  simplified.  Letting 
x  =  Q(r,r),  the  resulting  expression  is 


/C(r,r,Q(f,r))  =  ^lQ(r,Tj)  +  ^QT{r,r)Ufj 

-  -^(lQ{.r,f))(iQ{r,r))  -  3 (71) 

Thus, 

9fe(A:(r,  r,  x))  =  3?e{^/Q(r,  fj)  +  ^Q1  {r,r)Ufj}.  (72) 

The  following  formula  for  IC{r,f,Q(r,f))  is  obtained  in  a  similar  fashion.  The  fact  that  r  and 
f  are  eigenvectors  of  L-1  is  employed  in  the  computation. 


£(r,r,g(r,f))  =  ^3?e{/Q(r,£)  +  ^QT(r,  r)II£}- 
Condition  (H7)  may  now  be  rewritten  explicitly  as 


(73) 


&  +  A„(n)  <  0, 

where  /?2  is  as  defined  above,  and  where  A//  is  given  by  the  real  number 

A//(n)  :=  Q^rfUL-'fn  +  nL-^Qir,?) 


(74) 


1 


-  ~QJ  (r,  r){n(L  +  2iiocI)-1  +  ((L  -  2iu,cI)~1)TU}Q(r ,  f). 


(75) 
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It  is  not  difficult  to  ascertain  that  A//  is  nonnegative,  although  this  is  not  an  essential  consider¬ 
ation.  (One  proof  of  this  uses  the  fact  that  formula  (75)  is  a  special  case  of  the  sum  of  a  Hermitian 
form  and  a  quadratic  form.) 

We  have  just  rederived  the  following  known  criterion  for  asymptotic  stability  in  Case  (II). 


Theorem  3.  Let  hypothesis  (II)  hold,  and  suppose  that  fa  <  0,  where  fa  is  given  by  Eq.  (68). 
Then  the  origin  of  (1)  is  locally  asymptotically  stable. 

Regarding  the  structural  coefficients  and  that  have  not  been  specified  explicitly  in  the 
foregoing  analysis,  only  two  constraints  remain:  The  first,  the  symmetry  requirement,  entails  that 
the  value  of  a  coefficient  is  independent  of  the  order  of  subscript  indices.  r!  fie  second  constrain! 
is  that  the  function  V(z)  must  be  real-valued  for  x  £  IRn ■  Using  Proposition  3,  it  is  found  that 
this  latter  requirement  is  equivalent  to  what  might  be  called  a  conjugate  symmetry  relationship 
among  the  structural  coefficients,  the  exact  nature  of  which  is  specified  in  the  next  Corollary  to 
Proposition  3. 

Corollary  2.  (Conjugate  Symmetry  of  Structural  Coefficients  in  Case  (H))  Denote,  for 
any  positive  integer  i,  the  quantity  [i]  (the  “complement,  of  i ”) 


W  := 


i: 


if  i  =  1, 
if  i  =  2, 


(76) 


1  i  otherwise 

Then,  V(x)  is  real-valued  for  each  x  £  IRn  if  and  only  if  the  structural  coefficients  Kt^  and  Tt^:p 
satisfy  the  following  relationship: 


Kijk  ~  K[.'][jP]’ 

=  r[,][jP][P]-  (78) 

Proof.  Follows  immediately  from  Proposition  3.  c 

The  foregoing  construction  of  a  family  of  Liapunov  functions  V(z)  of  the  form  (59)  for  the  case 
in  which  D/( 0)  possesses  a  pair  of  pure  imaginary  eigenvalues  is  summarized  in  the  next  result  and 
algorithm. 


Theorem  4.  Let  hypothesis  (H)  hold,  and  suppose  that  fa  <  0,  where  fa  is  given  by  Eq.  (68). 
Then  any  function  V(x)  resulting  from  Algorithm  V//  below  is  a  Liapunov  function  for  the  equilib¬ 
rium  point  0  of  Eq.  (1). 
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Algorithm  V//.  (Construction  of  Liapunov  functions  V(x)  =  xTVx  +  K{x,  x,  x)  +  T (x,  x,  x,  x)  in 

Case  (II)) 

Step  1.  Compute  /  and  r  and  choose  a  basis  {ra,r2, . .  .  ,  rn)  for  IRn  with  r1  :=  r,  r2  :=  r,  and  for 
which  {r3, . .  . ,  r"}  is  a  basis  for  Es.  Calculate  the  row  vectors  {/3,  . . . ,  /n}  of  the  associated 
dual  basis.  Check  that  /?2  <  0,  where  is  as  defined  in  (68). 

Stop  2.  Pick  any  II  satisfying  (i)  Ilr  =  Ilr  =  0,  (ii)  inTIIir  >  0,  and  (iii)  wT(LTU  +  IlL)u)  <  0  for 
all  w  €  Es,w  ■£  0,  and  such  that  A//(II)  <  |/?2|,  where  A#(n)  is  given  by  Eq.  (75). 

Step  3.  Set  V  =  U  +  lTJ+ITl. 

Step  4.  Set  /C(x\x2,x3)  -  £?=]  E”=i  ELi  nijk(l'x')(P x2)(lkx3),  and 


Kill  = 

3^t'0(r'r)' 

(79) 

Kll2  = 

{4lQ(r,r)  +  2/Q(r,r)}, 

oiioc 

(80) 

Kill  = 

-\{4 lQ(r,  (2iucI  +  L)-V)  +  2 QT(r,  r)H(L  +  2iucI)~l r{},  i  =  3, . . 

• ! n! (81) 

K.12  = 

{lQ(f,  L- V)  +  lQ(r,  L- V)  +  QT{r ,  rjniTV'},  *  =  3, . . . ,  n. 

(82) 

Step  5.  Set  T(x\x2,x3,x4)  =  i  £?=i  £2=1  £p=l  ^(/‘x1  )(Px2)(/*x3)(Fx4).  Here,  rmi  (= 
T(r,  r,  r,  r))  and  run  (=  T(r,r,r,r))  are  selected  according  to  (H5)  and  (H6),  respectively. 

Step  6.  All  structural  coefficients  and  rtj*.p  which  have  not  been  specified  in  Steps  1-5  are 
assigned  arbitrarily,  modulo  the  symmetry  requirement  and  the  conjugate  symmetry  require¬ 
ment  (Eqs.  (77),  (78)). 
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8  CONCLUDING  REMARKS 


Liapunov  functions  for  nonlinear  systems  with  either  of  the  two  simplest  critical  cases  have  been 
explicitly  constructed.  For  the  case  in  which  the  sytem  linearization  possesses  a  simple  zero  eigen¬ 
value,  generically  the  Liapunov  functions  need  contain  only  quadratic  and  cubic  terms  in  the  state. 
However,  when  a  complex  conjugate  pair  of  simple,  pure  imaginary  eigenvalues  are  present,  the 
Liapunov  functions  contain  quartic  terms,  in  addition  to  the  quadratic  and  cubic  terms.  These 
Liapunov  functions  were  shown  to  predict  local  asymptotic  stability  precisely  when  certain  known 
sufficient  conditions  from  bifurcation  analysis  are  satisfied.  We  have  obtained  parametrized  “fam¬ 
ilies”  of  Liapunov  functions  for  the  studied  critical  cases,  in  the  same  sense  that  the  Liapunov 
matrix  equation  yields  an  infinite  set  of  quadratic  Liapunov  functions  for  asymptotically  stable  lin¬ 
ear  time-invariant  systems.  The  Liapunov  functions  are  computed  directly  in  terms  of  the  Taylor 
series  expansion  of  the  vector  field  f(x),  and  are  thus  amenable  to  symbolic  computer  coding.  The 
use  of  these  Liapunov  functions  in  the  design  of  feedback  control  laws  for  critical  nonlinear  systems 
is  a  topic  for  future  investigation. 

APPENDIX  A.  SOLUTION  OF  LINEAR  ALGEBRAIC  EQUATIONS 
WITH  SINGULAR  COEFFICIENT  MATRIX 

Consider  the  system  of  linear  equations 

Ax  =  b  (A.l) 

where  A  is  a  real  n  x  n  matrix  and  b  £  lRn.  Suppose  that  A  has  a  simple  zero  eigenvalue.  Let 
r  and  l  denote  right  (column)  and  left  (row)  eigenvectors  of  A,  respectively,  corresponding  to  the 
zero  eigenvalue,  and  require  that  these  be  chosen  to  satisfy  Ir  =  1.  Under  these  conditions,  the 
Fredholm  Alternative  asserts  that  (A.l)  has  a  solution  if  and  only  if  lb  =  0.  Moreover,  the  Fredholm 
Alternative  also  implies  that,  if  (A.l)  has  a  solution  x° ,  then  the  totality  of  solutions  is  given  by 
the  one-parameter  family  x  =  x°  +  ar  where  a  £  IR  is  arbitrary.  The  solution  is  rendered  unique 
upon  imposing  a  normalization  condition  which  specifies  the  value  of  lx. 

Introduce  subspaces  EC,ES  C  IRn  as  follows:  Ec  is  the  one-dimensional  subspace 

Ec  :=  spanjr},  (A. 2) 
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and  Es  is  the  (n  —  l)-dimensional  subspace 


EB  :=  {x  e  IRn\  lx  =  0}.  (A. 3) 

From  the  foregoing,  we  have  in  particular  that  if  lb  =  0  then  the  system  Ax  =  fc,  lx  —  0  has  a 
unique  solution.  Equivalently,  (A.l)  Aas  a  unique  solution  in  Es  for  any  vector  b  £  Es .  This  proves 
that,  the  restriction  [8,  p.  199]  A\E s  of  the  linear  map  A  to  Es  defines  an  invertible  (one-to-one  and 
onto)  map.  In  the  next  result,  we  exhibit  the  unique  solution  which  lies  in  Esoi  the  system  Ax  —  b, 
lx  —  0.  The  proof  is  elementary  [2]. 

Proposition  A.l  The  unique  solution  of  Ax  =  b,  lx  =  0  given  that  lb  =  0  is 

x  =  {ATA  +  lTiy'ATb.  (A. 4) 

This  result  motivates  the  following  introduction  of  notation: 

A~  :=  (ATA  +  lTl)-'lAT.  (A. 5) 

Thus,  the  inverse  of  the  restricted  map  A\E>  exists  and  is  given  by 

(A\E.)-l  =  A-.  (A.6) 
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